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The high T c cuprates have a complex phase diagram with many competing phases. We propose a 
bosonic effective quantum Hamiltonian based on the projected SO(5) model with extended interac- 
tions, which can be derived from the microscopic models of the cuprates. The global phase diagram 
of this model is obtained using mean-field theory and the Quantum Monte Carlo simulation, which is 
possible because of the absence of the minus sign problem. We show that this single quantum model 
can account for most salient features observed in the high T c cuprates, with different families of the 
cuprates attributed to different traces in the global phase diagram. Experimental consequences are 
discussed and new theoretical predictions are presented. 
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I. INTRODUCTION 

On the first look, the phase diagram of the high transi- 
tion temperature superconducting (HTSC) cuprates has 
a striking simplicity: there are only three universal 
phases in the phase diagram of all HTSC cuprates: the 
antiferromagnetic (AF), the superconducting (SC) and 
the metallic phases, all with homogeneous charge distri- 
butions. However, closer inspection shows a bewilder- 
ing complexity of other possible phases, which may or 
may not be universally present in all HTSC cuprates. A 
large class of these phases have inhomogeneous charge 
distributions. Because of this complexity, formulating 
an universal theory of HTSC is a great challenge. The 
50(5) theory unifies the AF and the SC order parame- 
ters into a single five dimensional order parameter called 
the superspin, and the effective quantum theory of the 
superspin naturally explains proximity between the AF 
and the SC phases in the observed phase diagram 1 .. The 
Goldstone modes of the superspin fluctuations can be 
identified with the it resonance mode observed in the neu- 
tron scattering experimento 2 i 3 i 4 i 5 i 6 i 7 i 8 i 9 i 10 i 1:L i 12 i 13 i 14 i 15 i 16 . 
This theory also predicts the AF vortex stat e 1 ' 17 , 
which has recently been observed in a number of 
experimenta 18 ! 19 . 20 ' 21 . 22 ! 23 . 24 ! 25 ' 26 . Initially, the SO(5) 
theory was motivated by the simplicity of the pure AF 
and SC states, however, given the encouraging agree- 
ments with the experiments, it is tempting to construct 
an unified theory of the global phase diagram of the 
HTSC which addresses the more complex inhomogeneous 
phases as well. Complexities can of course be introduced 



phenomenologically into the Landau-Ginzburg type of 
theories by simply introducing more order parameters. 
However, this type of approach necessarily limits the pre- 
dictive power of theory. The goal of this paper is to 
present a single effective quantum model of the superspin 
degree of freedom, which can be derived systematically 
from the microscopic electron models, and can be inves- 
tigated reliably both analytically and numerically. The 
global phase diagram of this model is then compared with 
the experimentally observed phase diagram of the HTSC 
cuprates. 

When formulated on a coarse-grained lattice, with high 
energy charge states projected out, the projected SO (5) 
model describes five local superspin degrees of freedom 
per plaquetteS 7 .. These five states are the spin singlet 
state at half-filling, the spin triplet states at half-filling, 
and the singlet d-wave hole pair state. Using the Con- 
tractor Renormalization Group (CORE) algorithm, Alt- 
man and AuerbachS showed that the projected 5*0(5) 
model can be systematically derived from the microscopic 
electron models, and they also determined the parame- 
ters of the effective SO (5) model explicitly from the mi- 
croscopic interaction parameters (see also Ref. 29]). Re- 
stricted within the subspace of these five local states, 
the Hamiltonian describing their propagation and inter- 
action is completely expressed in terms of bosonic oper- 
ators and can be studied reliably by the quantum Monte 
Carlo (QMC) calculations. The simplest form of the pro- 
jected 5*0(5) model has been studied extensively by the 
QMC method both in two dimensions^iSii 3 ^ and in three 
dimensions 33 . The overall topology of the phase diagram, 
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the scaling properties near the multi-critical point and 
the nature of the collective excitations can be reliably 
obtained from the QMC method, within the parameter 
regime of experimental interests. 

The simplest form of the quantum 50(5) model de- 
scribe either the direct, first order transition from the AF 
to the SC state, or two second order transitions with an 
uniform, intermediate AF/SC mix phase in between the 
pure AF and the SC statesiSl. In the case of the direct 
first order transition as a function of the chemical po- 
tential, the system at a fixed density is phase separated. 
However, in the HTSC cuprates, there are other forms 
of charge and spin ordered states. For example, neu- 
tron scattering cross section in LSCO material is peaked 
around (tt ± 5, tt) and (tt, tt ± S), where S ~ l/8 34 i 3 ^ 3 9. 
STM experiments have revealed periodic charge modula- 
tion with period close to four lattice spacingi^S, either 
near the vortex core or near surface impurities. In lat- 
ter case, alternative interpretatio n 39 ! 40 ' 41 based on quasi- 
particle interference is also possible, and the two points 
of views are summarized by Kivelson et al.^ 2 . Motivated 
by these experiments, we extend the simplest form of the 
projected 5*0(5) model to include extended interactions 
among the five bosonic states. In fact, these extended in- 
teractions also arise naturally by carrying out the CORE 
algorithm to extended ranges. 

The projected 50(5) model with extended interac- 
tions support a more complex phase diagram. In par- 
ticular, there are insulating phases at fractional filling 
factors where the charges form a lattice, usually com- 
mensurate with the underlying lattice. A crucial aspect 
of this model is that all charge density wave states are 
formed by the Cooper pairs of the holes, rather than the 
holes themselves 4 ^. Throughout this paper, we shall de- 
note such states as the pair-density-wave (PDW) states 
or pair checkerboard states. This distinction has a pro- 
found experimental consequence, since the real space pe- 
riodicity of the former is larger than the latter by a factor 
of y/2. This type of insulating PDW states is a conse- 
quence of strong pairing and low superfluid density, a 
condition which is naturally fulfilled in the underdoped 
cuprates, but has not yet been unambiguously identi- 
fied in other experimental systems before. The PDW 
state can either take the form of stripes or checkerboards, 
depending on the ratios of the extended interaction pa- 
rameters in the model. Furthermore, PDW states with 
longer periodicity generally requires longer range inter- 
actions to stabilize. Based on this reasoning, a simple 
picture emerges for the global phase diagram of under- 
doped cuprates. The phase diagram consists of islands 
of insulating PDW states, each with a preferred rational 
filling fraction, immersed in the background of SC states 
(see Fig. 0J. The height of the Mott insulating PDW 
lobes vary depending on the preferred filling fraction and 
the range of extended interactions, but in principle, these 
insulating states are all self-similar to each other, and 
similar to the parent AF insulator at half-filling. There 
can be either a direct first order phase transition, or two 



second order phase transitions between the SC state and 
the PDW state, with the possibility of an intermediate 
"supersolid" phase, where both orders are present. 

Based on our model, the bewildering complexity of the 
cuprate phase diagram can be deduced from a simple 
principle of the "Law of Corresponding States". This 
concept is borrowed from the work of Kivelson, Lee and 
Zhang on the global phase diagram of the quantum Hall 
effect 4 ^, in fact, our proposed phase diagram in Fig. 0] 
bears great similarity to Fig. 1 of that reference. In 
the case of the QHE, the "Law of Corresponding States" 
physically relates all quantum phase transitions at vari- 
ous filling fractions to a single quantum phase transition 
from the v = 1 integer state to the Hall insulator. In 
recent years, this powerful mapping among the differ- 
ent fractional states has been made more precise by the 
derivation of the 5L(2, Z) discrete modular group trans- 
formation from the Chern-Simons theorji 4 ^ 4 ^ 4 ^. Simi- 
larly, the central idea of the current paper is to relate 
the fractional Mott insulator to SC transition with the 
transition from the AF Mott insulator at half-filling to 
SC state, which is already well understood within the 
context of the original, simple 50(5) theory. The con- 
struction of the Mott insulating states at various frac- 
tional filling factors can be constructed from the "Law 
of Corresponding States" , iterated ad infinitum, to give 
a beautiful fractal structure of self-similar phases and 
phase transitions, as presented in Fig. 0] The various 
different compounds of the HTSC cuprates families have 
slightly different microscopic parameters, and they cor- 
respond to different slices of this global phase diagram. 
The global phase diagram provides a basic road-map to 
understand the common elements and differences among 
various HTSC compounds. 

This paper mainly focuses on the zero temperature 
global phase diagram of the underdoped cuprates. How- 
ever, it is understood that the model is valid below the 
pseudogap temperature, which we interpret as the tem- 
perature below which the system can be effectively de- 
scribed by the collective bosonic degrees of freedom, like 
the magnons and the hole pairs. Therefore, it is implied 
that the pseudogap state is a regime where the various 
ground states discussed here compete with each other, 
and different experiments may access different aspects of 
these competing states. The existence of the pseudogap 
temperature gives the fundamental experimental justi- 
fication to investigate the global phase diagram of the 
underdoped cuprates by a purely bosonic model. In the 
future, we shall use the same model to investigate the 
manifestations of these competing states at finite tem- 
perature, in the pseudogap regime. A comparison of the 
charge order predicted by this work and the STM experi- 
ment in the pseudogap regime has recently been reported 
in Ref.|H. 

While this paper is presented within the logical context 
of the 50(5) theory, some of the ideas and results bear 
intellectual similarities to the previous theoretical works. 
The idea of doped holes forming ordered stripes has been 



discussed extensively in Ref. |34l49l50l5ll52l53j . Al- 
though we focus more on the charged ordered states in 
the forms of checkerboards of hole pairs, they are concep- 
tually related to stripes and can be realized experimen- 
tally or theoretically depending on the microscopic pa- 
rameters. The pseudogap temperature was identified as 
the formation temperature of Cooper pairs by Emery and 
Kivelson^a. Our interpretation of the pseudogap temper- 
ature is more general, which also includes the formation 
of the magnetic collective modes in addition to the holes 
pairs. Vojta and Sachde^iS have discussed the phase 
diagram of doped Mott insulator with various charge or- 
dered insulating states at rational fractions. More re- 
cently, Zhang, Demler and Sachdev have studied exten- 
sively the competition among charge and spin ordes££*£t 
Laughlin pointed out that the small superfluid density in 
the underdoped regime is responsible for various charge 
ordering phenomena—. Haas et alJ& have noticed that 
the Wigner crystal state of the hole pairs could be sta- 
bilized due to the competition of phase separation and 
long ranged Coulomb interaction. Kim and Hor 60 have 
discussed experiments at certain "magic" filling fractions 
in terms of the commensurate Wigner crystal type of or- 
der of the electrons, rather than the hole pairs discussed 
in this paper. Restricted to the charge sector, the pro- 
jected £0(5) model is essentially the same as the hard- 
core quantum boson model on a lattice, whose phase di- 
agram has been extensively studieo&SSiS. 

This paper is organized as follows. In section [HI the 
projected 50(5) model with extended interactions is pre- 
sented. The choice of parameters is discussed from the 
CORE algorithm and phenomenology. The self-similarity 
of the insulating states and the classification of the quan- 
tum phase transitions are then discussed in section ITTfl 
In section IIVI the global phase diagram of the model 
is obtained within the mean field theory. The low en- 
ergy collective modes and their quantum symmetry is 
then studied using a slave-boson approach. In section Ivl 
QMC simulation is carried out to compare with mean- 
field results obtained in section IIVI The experimental 
consequences and predictions are discussed in section lvTl 
Finally, section IVTT1 concludes our study. 
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II. HAMILTONIAN OF THE MODEL 
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FIG. 1: Illustration of the basic idea of the CORE method. 
To implement the CORE method, first decompose the orig- 
inal lattice in plaquettes, and then truncate the spectrum 
of a given plaquette to five lowest states, i.e, the singlet, 
hole-pair and three magnon states. An effective Hamiltonian 
for these bosons can then be determined using the CORE 
method. Left: local bosons in the original lattice. Gray rect- 
angle denotes the singlet RVB vacua, circles denote holes and 
the set of two parallel vertical arrows denote the magnon. 
Right: local bosons on the lattice of plaquette. Leaf-like pat- 
tern denotes a local d-wave hole-pair on a plaquette. Canted 
arrow denotes local magnon on a plaquette. The singlet RVB 
vacuum is denoted by an empty site. 



spectively. The low energy eigenstates of the Heisenberg 
plaquette can be determined easily. The nondegenerate 
ground state |fi) (see Ref.[28| for a real-space represen- 
tation in terms of the microscopic states on a plaquette) 
has energy E = -2J and total spin 5 = 0. This "RVB" 
like singlet state will be the vacuum state of the effective 
bosonic 50(5) model. The next set of energy eigenstates 
are three triplet states tjjfi) with energy E t — ~J and to- 
tal spin quantum number 5=1. All other energy eigen- 
states of the Heisenberg plaquette have energies E > 
and can be neglected in the low energy effective model. 
It should be noted that the operator P a with spin 1 and 
charge create hardcore bosons because one cannot cre- 
ate more than one of them simultaneously on a single 
plaquette. The ground state of two holes is a "Cooper" - 
like hole pair with internal d-wave symmetry with respect 
to the vacuum. 

Using the CORE method and keeping only the five 
lowest states (the singlet boson, the three magnons t\ 
and the hole-pair t\), the effective Hamiltonian of these 
bosons can be obtained a£ 28 - 29 - 64 



The effective bosonic SO (5) model can be derived 
directly from the microscopic Hubbard model or t- 
J model, through a renormalization group transfor- 
mation called the Contractor Renormalization(CORE) 
method22iS4. To construct bosonic quasiparticles from 
fermionic model, we divide the lattice into effective sites 
containing an even number of sites. In order to conserve 
the symmetry between x and y-direction in the system, 
a plaquette of 2 x 2 sites are typically chosenSLS. In 
their CORE study of the 2D Hubbard model, Altman and 
Auerbach2& started from the spectrum of lowest-energy 
eigenstates of the 2x2 plaquette for 0, 1 and 2 holes, re- 



Ji. = 7io + Tie 



(1) 



where 7io is the Hamiltonian of the previously studied 
50(5) model containing only on-site interaction s 27,30 ! 33 

Ho = A c E4(0*h(i) + A. EE 

i a i 

- JoJ2[tl(i)t h (j)+H.c] 

- EE [*«(*)+*«(*)] +*«(?)]. (2a) 

a (ij) 
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and TLext is the part containing extended interactions 



He 



nh(i)n h {j) 

- J -££ [4(0*/»Cj)4(j)*aW+-ff.' 

+ ^r££["a(i)«a(i) + n/ v ^')n a (i)] 

+ E E ^s(^)s(^)s- (2b) 

(ij) 5=0,1,2 

The model is subjected to the hard-core constraint 
£n a (i) < 1. 



(3) 



Here, A c and A s are the energy costs to create a hole- 
pair and magnon respectively. J c and J s are the hopping 
terms of hole-pairs and magnons. th{i) and tl(i) are 
the annihilation and creation operators of hole-pair on 
plaquette i. t a (i) and v a (i) are the annihilation and cre- 
ation operators of magnon on plaquette i for a = x, y, z. 
(■ ■ • ) and ((■ ■ • )) denote nearest- neighbor (nn) and next- 
nearest- neighbor (nnn), respectively. rih(i) = ^ 
and n a (i) = t\{i)t a {i) are the hole-pair density and 
magnon density operators on plaquette i, respectively. 
The hole-pair density of nh per plaquette corresponds to 
twice the real doping of S holes per lattice site, i.e. 



S = n h /2. 



(4) 



Finally, (t»ij)o creates two magnons simultaneously on 
plaquettes i and j, which are coupled into total spin 
S. Ti-ext contains nn and nnn hole pair interactions V c 
and V£, exchange hopping terms J v between t a and th 
bosons, interaction between nn magnons and hole- 
pairs. The 4-magnon interactions Vb,i,2 are important 
in the pure AF phase but we can neglect them since we 
are mostly interested in the doped phase where magnon 
density decreases. According to the CORE calculation 
on 2-plaquettes and fixing J s = 1 as the unit of energy, 
one obtains from the t-J model with J/t ~ 0.4 (relevant 
for cuprates) : 
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For J c ~ 2, H.q is approximately SO (5) symmetric at the 
mean field level 27 i 30 i 33 . 

The CORE derivation of the 5*0(5) model Q is only 
approximate. It should also be born in mind that the t—J 
and the Hubbard models are also approximate models 
of the real cuprates themselves. If one started from a 
different microscopic model (for example, next-nearest 
neighbor hopping, extended Coulomb repulsion and etc), 
one would have obtained a similar effective Hamiltonian 



with different parameters. Therefore, in this paper, we 
shall take the CORE parameters as a guidance, and study 
the robust properties of the SO (5) model with a more 
general set of parameters, as to reproduce well-known 
results and compare directly with experiments. 

At half-filling [nu — 0), the model involves only the 
singlet and the magnon, and the effective Hamiltonian 
containing only J s and A s can be rewritten as 

rr _ o j \^ 4 (») +*<.(*) tUj)+t*(j) A s ^ 2 

H - 2J S ^ -j= + — 2^ L a fsV 1 )^ 



V2 



(5) 

where [tj,(«) +t a (i)}/^/2 is the AF moment and L a p(i) is 
the 50(3) symmetry generato»21 



L al3 (i) = -iltt(i)te(i)-tl(i)t a (i)] 



(6) 



This model is similar to the nonlinear a model 
(NLcrM)&£& 



(7) 



where mf is the a component of the AF moment and Si 
is the angular moment on site i. After rescaling of time 
using the spin velocity c = y/PsJx an d up to a prefactor, 
the Lagrangian density of the NLtrM can be cast in the 
usual form in the continuum limit 



Cnlct = -(d x m) 2 + g(d T m) 2 , 
9 



(8) 



where g ~ 1 1 J p s X is a dimensionless constant. This 
model has been studied extensively^. It has a transition 
towards a disordered state at g c = 1.45. From the com- 
putation of the staggered moment, we can find g such 
that the original Heisenberg value (0.3) of the AF mo- 
ment is recovered: gn = 1.125. On the other hand, we 
know from mean-field calculations and QMC simulation 
that the disordered phase occurs at A s /J s = 8. We then 
obtain the proportionality factor between g and 1/ ' -JpsX- 
Using gn , we find that an effective model for Heisenberg 
corresponds to A s /J s = 4.8. 

In most parts of this paper, we shall consider the sim- 
plified model with = V w = 0. 



III. HEURISTIC ARGUMENT ON THE 

SELF-SIMILARITY BETWEEN 
PAIR-DENSITY- WAVE (PDW) STATES 

Let us ignore for a moment the magnons, and consider 
a hard-core boson model with extended interactions. The 
phase diagram£i*£2 of a hard-core boson model with nn 
interaction contains one superfluid state and three Mott- 
insulating states, corresponding to zero-doping (n^ = 0), 
half- filling (n^ = 1/2) and fully-occupied (rih — 1), as 
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FIG. 2: Phase diagram of a hard-core boson model with 
nearest-neighbor interaction. There are one superfluid (SF) 
state and three insulating states: zero-doping state (rih ~ 0), 
half-filled state (nj = 1/2) and fully-occupied state (rih = 1). 
At the next level of the hierarchy, longer ranged interactions 
lead to new insulating states with nj, = 1/4 and nu = 3/4, as 
shown in Fig. |3] 

sketched in Fig. The half-filling state has checker- 
board charge order. If one transforms the hard-core bo- 
son model into an AF Heisenberg model, the checker- 
board order of the bosons simply correspond to the AF 
order of the Heisenberg spins. The following argument 
assumes that the checkerboard order at half-filling is a 
basic and robust form of order, such that it is repeated 
at all different levels of the hierarchy. 




FIG. 3: Phase diagram of a hard-core boson model with nn 
and nnn interactions. There are a superfluid (SF) state and 
five insulating states with doping 0/4, 1/4, 2/4, 3/4 and 4/4. 
At the next level of the hierarchy, new insulating states are 
developed at doping level nj, = 1/8, = 3/8, = 5/8 and 
rih = 7/8. This hierarchy construction can be iterated ad in- 
finitum, to obtain a self-similar phase diagram with insulating 
phases at doping level p/2 n , with integers p and n, such that 
0<p<2". 

If we regard the empty sites of the half-filled checker- 
board state as an inert background, we obtain a fully- 
occupied Mott-insulating state on the coarse-grained lat- 
tice, with lattice spacing y/2a x s/2a. The nnn interac- 
tion on the original lattice becomes the nn interaction on 
the coarse-grained lattice, and a new half-filled checker- 



board state can be stabilized on the coarse-grained lat- 
tice. Such a state corresponds to doping rih = 1/4 on 
the original lattice. Similarly, we can regard the filled 
sites of the original half-filled checkerboard lattice as an 
inert background, leaving with an empty state on the 
coarse-grained lattice. A new checkerboard state can 
again be stabilized on the coarse-grained lattice, which 
corresponds to doping rih = 3/4 on the original lattice. 
This hierarchical procedure of forming a new daugh- 
ter checkerboard state from two parent checkerboard 
state can obviously be iterated ad infinitum, to obtain a 
fractal-like, self-similar phase diagram as shown in Fig|3] 
It is interesting to note that the nearest-neighbor interac- 
tion on the coarse-grained lattice is just the next-nearest- 
neighbor interaction on the original lattice. There is 
also the possibility that small regions with coexisting SF 
and PDW orders ( "supersolids" ) are present around the 
Mott-insulating lobes in the phase diagramSiSiSLS. 

Having presented the generic phase diagram for the 
charge boson only, we consider now the inclusion of the 
magnons in the full SO (5) with extended interactions. 
Generally, charge ordered insulating states also have AF 
order. The rih = state of the charge boson corresponds 
to the undoped parent Mott insulator. The rih = 1/2 
state of the charge boson would correspond to S = 1/4 
doping for the cuprates, which is probably at or beyond 
the limit of applicability of our bosonic model. There- 
fore, the phase diagram of the hard-core boson model in 
the range of < rih < 1/2 from Fig. [3] would translate 
into a phase diagram of the cuprates in the doping range 
of < S < 1/4, as shown in Fig. 0] As we shall show 
later, this phase diagram is supported by accurate QMC 
calculations of the SO(5) model with extended interac- 
tions. We expect that the insulating states of 8 = 1/16 
and 1/8 are AF ordered. Since the magnon density de- 
creases with increasing doping, the S = 3/16 state may 
not be AF ordered. 

The nature of the phase boundary between two differ- 
ent phases shown in Fig. 01 requires careful characteriza- 
tion. We can classify all phase transitions into two broad 
classes. Class A describes transitions at fixed chemical 
potential, typically at an effectively particle-hole sym- 
metric point around the tip of the Mott lobe. Class B 
describes transitions where the chemical potential or the 
density is varied. Each broad class is further classified 
into three types, 1, 2 and 1.5. Generically, the phase 
transition between two ordered phase can be either a 
single first-order transition or two second order transi- 
tions, with a mixed state in between, where both order 
parameters are nonzero. A third marginal possibility oc- 
curs at a symmetric point, when these two second order 
phase transitions collapse into a single one. In the con- 
text of high Tc cuprates, these three types are shown 
in phase diagrams of Ref.pJ as Fig. la,b,c, respectively. 
This situation can be easily understood by describing the 
competition in terms of a Landau-Ginzburg functional of 
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FIG. 4: A typical global phase diagram of the cuprates 
in the parameter space of chemical potential and the ratio 
of boson hopping energy over Coulomb interaction energy. 
This phase diagram shows self-similarity among the insulating 
PDW states at half-filling and other rational filling fractions. 
There are two types of superfluid-insulator transition. The 
quantum phase transition of "class A" can be approached by 
varying the hopping energy, for example, by applying a pres- 
sure and magnetic field at constant doping. The quantum 
phase transition of "class B" can be realized by changing the 
chemical potential or doping. There can be either a direct first 
order phase transition, or two second order phase transitions 
between the SC state and the PDW state, with the possibility 
of an intermediate "supersolid" phase, where both orders are 
present. Different families of cuprates correspond to different 
traces of "class B". For example, we believe YBCO is Bl- 
like, BSCO may be close to B2-like and LSCO is B3-like. 
The vertical dash-dot line denotes a boundary in the over- 
doped region beyond which our pure bosonic model becomes 
less accurate. 



two competing order parameters^, which is given by 



F 



-ri4>i + + U X 4>\ + U 2 (f>t + 2Ml20!</»2 (9) 



where <pi and <p2 are vector order parameters with iVi and 
N2 components, respectively. In the context of 5*0(5) 
theory, N\ = 2 and N2 = 3, and we can view <j)\ as SC 
component of the superspin vector, and as the AF 
component of the superspin vector. These order param- 
eters are obtained by minimizing the free energy F. By 
tuning t 1 , one can drive a quantum phase transition from 
AF to SC. For U12 > ^/uiU2, the quantum phase tran- 
sition from AF to SC is a single first order transition of 
"type 1" . For U12 < ^/ 1*11*2, the transition from AF to 
SC consists in two second order transitions, and there is 
a finite range of r\ where AF and SC coexist uniformly; 
the transition is of "type 2" . For u\2 — ^Ju\U2, the phase 
transition occurs at 



where the free energy takes the 50(5) symmetric form 

(11) 



TP r W?!± 12 , 2 74 

F = — - — 4> +u 12 4> 



with 



72 _ 4>l <t>2 



Ui *JU2 



(12) 



'•1 



''2 



(10) 



Since the free energy depends only on </>, one order pa- 
rameter can be smoothly rotated into the other without 
any energy cost. At this point, the chemical potential is 
held fixed, but the SC order parameter and the charge 
density can change continuously according the condition 
that 4> 2 is constant. This is a special case intermediate 
between "type 1 and 2" transitions, where two second 
order phase transitions collapse into one. This transition 
can only occur at a 5*0(5) symmetric point. We thus 
classify it as "type 1.5". The full quantum 50(5) sym- 
metry can only be realized in the class A transition of 
"type 1.5". On the other hand, the static, or projected 
50(5) symmetry can be realized in class B transitions of 
"type 1.5". 

In HTSC cuprates, the charge gap at half-filling is very 
large, of the order of U ~ 6eV , it is not possible to induce 
the "class Al" transition from the AF to the SC state 
by conventional means. However, the charge gap in the 
fractional insulating states is much smaller, of the order 
of J c , and it is possible to induce the "class A2" or "class 
A3" insulator to superconductor transition by applying 
pressur e) 70 ! 71 or by applying a magnetic fielcS 2 * 7 ^. 

As the chemical potential or the doping level is var- 
ied, a given system, roughly corresponding to a fixed 
value of the quantum parameter J c /V c , traces out dif- 
ferent one dimensional slices in this phase diagram, with 
typical slices Bl, B2 and B3 depicted in Fig. The 
nature of the phase transition Bl is similar to that of 
the superspin- flop transition discussed in Ref.[lJ. In this 
case, the phase transition from the AF to the SC state 
can be further classified into "types 1, 1.5 and 2", with 
the last two cases leading to a AF/SC mixed phase at 
the phase transition boundary. For lower values of J c /V c 
, the trace B3 encounters the S = 1/8 insulating phase. 
The key signature of this type of phase transition is that 
the SC T c will display a pronounced minimum as the dop- 
ing variation traces through the (5=1/8 insulating state. 
Meanwhile, the AF ordering (possibly at a wave vector 
shifted from (7r,7r)) will show reentrant behavior as dop- 
ing is varied. The phase transition around the fractional 
insulating phases can again be classified into types "1, 
1.5 and 2", with possible AF/SC, AF/PDW, SC/PDW 
and AF/PDW/SC mixed phases. 

We believe that the AF to SC transition in the YBCO, 
BCCO and the NCCO systems corresponds to a "class 
Bl" transition. These systems only have a AF to SC 
transition, which can be further classified into types "1, 
1.5 and 2", but they do not encounter additional stati- 
cally ordered fractional insulating phases. On the other 
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hand, the phase transition in the LSCO system, where 
T c displays a pronounced dip at S — 1/8, correspond to 
the "class B3" transition. 



FIG. 5: The schematic plot of the unit cell of a quarter-filled 
magnetic insulating state. The solid square denotes a plaque- 
tte of the original lattice and the dash-dot square denotes a 
2-plaquette by 2-plaquette unit cell. 



IV. MEAN-FIELD PHASE DIAGRAM OF THE 
MODEL 

A. Four-sublattice ansatz and mean-field phase 
diagram 



Since the Hamiltonian Q contains up to nnn inter- 
actions, we can introduce the following four-sublattice 
ansatz within mean held theory: 



[E Ei 

i.B..Hi 

[EE] 
IE El 



IE E : 
|HEj 

[EE] 
IE El 



I*) 



+ h m t\{m,r) + x m tl(m,r) (13) 



where |fi) is the singlet ground state, e m , x n 
variational parameters, m — A, B,C, D denote the sites 
in a unit cell, and r is the coordinate in the lattice of unit 
cells, as sketched in Fig. The mean-field energy Emf 
reads 



N = A c (h\ + h%+hc + h 2 B ) - 2J C [e A h A + e c h c ] [e B /i B + e D h D ] + A s (x\ + x% + x 2 c + x 2 B ) 

-4J S [eAXA + e c x c ] \s-bxb + e D x D ] - 4J W [h A x A + h c x c ] + h D x D ] + 2V c [h 2 A + h 2 c ][h 2 B + h 2 D ] 

+Wt[h\h 2 c + h 2 B h 2 D ] + 2V^ {[h\ + h 2 c ][x% + x%] + [h\ + h 2 D \[x\ + x 2 c ]) (14) 



with the hard-core constraint 

Kn + x m + e li = !> m = A B, C, D. (15) 

Here, N/A is the number of unit-cells and N is the num- 
ber of plaquettes. 

By minimizing the energy functional l|14fl subjected 
to the hard-core constraint (|15|l . we obtain the mean- 
field ground state for a given set of parameters. Fig. 
plots the mean-field global phase diagram, for A s = 4.8, 
V c = 4.1010, VI = 3.6329 and J„ = V„ = 0. 

This phase diagram displays some rich features as ex- 
pected. It has three insulating state: an undoped anti- 
ferromagnetic (AF) state, an insulating AF PDW state 
with hole-pair density = 1/4 (5 = 1/8) and an in- 
sulating PDW state with hole-pair density = 1/2 
(5 = 1/4). Besides these insulating states, it also has 
a pure SC phase, a supersolid phase and mix phases of 
coexisting AF and SC order. 

In Fig. and Fig. |SJ we plot the doping dependence 
of SC and AF orders for different J c . If one follows a 
"class Bl" trace, such as the one with fixed J c = 1.5, 
the doping dependence of SC order mimics the behav- 



ior of YBCO and BSCO families with a underdoped 
region rih < 0.3(5 < 0.15) and an overdoped region 
rih > 0.3(6 > 0.15). If one follows a "class B3" trace, 
such as the one with fixed J c = 1.0, the SC order dis- 
plays a pronounce dip and the AF ordering is strongly 
enhanced around rih — 1/4 (6 — 1/8). Therefore, the 
"class B2" trace mimics the behavior of LSCO family. 

The doping dependence of charge order parameter is 
also plotted in Fig. |5J It measures the charge modulation 
defined by 



PDW =-Y^\n h {m)-n h \ (16) 



where is the average hole-pair density and m is 
summed over A, B, C, D. While "class Bl" trace shows 
no charge ordering in underdoped region, "class B3" 
trace displays a clear signature of charge ordering around 
5 = 1/8. 
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FIG. 6: The MF phase diagram obtained by minimizing the 
energy functional l|14[l subjected to the hard-core constraint 
1151 . J tv and Vn are taken to be zero. 



FIG. 8: Doping dependence of AF order parameter for differ- 
ent J c . For J c > 1.3, AF order decreases as doping increases 
and vanishes around 5 = 0.1. 




FIG. 7: Doping dependence of SC order parameter for differ- 
ent J c . For small J c , there is a dip around hole-pair doping 
rih = 1/4 (real doping 8 = 1/8). 



FIG. 9: Doping dependence of PDW order defined for differ- 
ence J c . For J c < 1.3, a strong peak is present at rih = 1/4 
or real doping 5 — 1/8. 



B. Slave-boson approach, effective Hamiltonian 
and dynamical 5*0(5) symmetry 

The hard-core constraint J3J) can also be enforced by 
introducing a slave-boson t\(i) for each latticed. The 
presence of this boson indicates that the plaquette i is 
empty. The hard-core condition is then replaced by 

$W = f E + + t%(i) - 1 J = 0. (17) 



This constraint can be enforced by introducing an addi- 
tional site dependent field A(i) and adding to the Hamil- 
tonian an additional term 

Hx = -J2WQ(i)- (18) 

i 

Since in physical states one always has one and only one 
boson per lattice sites, destruction (creation) of a boson 
t a (t' a ), a = x,y,z,h, must always be accompanied by 
creation (destruction) of the empty boson t\{t e ). In this 
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way, the whole Hamiltonian takes the form 
H sb = A c 5^tj;(i)t/i(i) + A.534(i)t a (i) 

i i 
-Jc52[A(i)te(i)ti(3)th(j)+H. 



-J s J2[tUe(i)+ttt a (i)] [tUe(j)+ttt a (j)] 



^A(z)Q(i)+W e 



(19) 



By integrating out the A field in the partition func- 
tion, one automatically enforces the hard-core constraint 
Q(i) — on each site. The saddle point approximation to 
(|19|l corresponds to replacing the boson operators and A 
field with real constants to minimize the energy Within 
the four-sublattice ansatz, this is equivalent to minimiz- 
ing the energy functional subjected to the hard-core con- 
straint (JTSJ. After obtaining the ground state, we can ex- 



pand the boson operators around their mean-field values 
and drop quartic terms in the Hamiltonian i|19f) to yield 
an effective Hamiltonian of the boson operators around 
the ground state. 

We shall now study the checkerboard state character- 
ized by following mean-field solution: 



i-c = no 







(20) 



for the simplified model with J w = V w = 0. 

The saddle-point of the slave-boson Hamiltonian 119|) 
can be solved to yield 



xd = x B = x c = 

gd = e B = e c = 1 

A_d = = Xb = Ac = 

X A = A c 



(21a) 
(21b) 
(21c) 
(21d) 



for A s /Js> and 



XD 



xc 



XB 



\ 



1 

Ic 



/64J| + A 2 
16 J? + A 2 



\ 



'16J S 2 + A 2 
64J? + A 2 



A. 64J? + A| 

2 ZJs \l6J 2 + A 2 



Xb = - 4 J s 



1 16 J 2 
64J? 



en = e_e = 



\ 



1 

Tg 



ec 



\ 



A.4 = A c 




64 J| + A 2 S 
16 J 2 + A? 



1 




/l6Jf- 


\- A 2 


8 




/ 64 J s 2 - 


1" A 2 _ 



(22a) 

(22b) 
(22c) 
(22d) 



for A s /J s < 4^/2. 

The bosons are then expanded around their mean-field 
values as 



t h (m 


r) 


= + 




(23a) 


t x [m 


r) 


— 


b x (m,r), 


(23b) 


t e (m, 


r) 


e m -f- 


b e (m,r), 


(23c) 


t y (m, 


r) 


= & y (m 


r), 


(23d) 


t z (m, 


r) 


= M m 


r). 


(23e) 



Again, m = A, B,C, D and r is the coordinate in the 



coarse-grained lattice of unit cells. Plug l|23a|) - (|23e|) into 
the slave-boson Hamiltonian (|19fl and drop the quartic 
terms to yield the quadratic effective Hamiltonian 

H eff = E Q + Hf f + H* ff + H e y !i + H e z ff (24) 
with Eq the mean-field ground state energy given by 

E a = 2N{A c -X B )x%+N{A c -Xc)x c 

+ANx c ecXBeB (25a) 

and //;. / / .//: // .//; f / and H^f 
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H 



eff 



H 



m=B,C,D 



E 

-B,C 

J *E 

q 

2 ^E 



^(A s - X m )b x b x {m, 9) - E ^mbthim, q) 
q q 

ie B x B (bl(C, q)b e (C, q) + H.c\ + B^C 
ec[bUC, q) + k(C, -q)} + x c [bt(C, q) + k{C, -q)] 



cos q x (e B [bUB,-q) + b x (B,q)]+XB[bl(B,-q)+b e (B,q)] + cos q y ( B <-> D 



(25b) 



eff 



E - X b+ 2V c) [bik(B, q) + b{b h (D, q)\ + £ (A c - A c + 4V C ') Sj^C, g) 
q q 

^(A s - A c )St b x (A, q)-X A J2 &M A > «) - J s E 8esa;s ^)&e(A g) + S|(A, g)S x (A g)] 



2J c e B 



cosq y ( S e (A, q)b h (B, -q) + b[(B, -q)b\{A, q) ) + cosqj B D 



2J c esec'^2 



cosq x b\{C, q)b h (B, q) + b\{B , q)b h (C, q) )+ cosgJ B D 



(25c) 



Hi" = ]T(A S - X B ) (bib a (B, q) + blb a (D, q)) + ^(A s - Xc)P a b a (C, q) + ^(A s - A c )P a b a (A q) 



-2J s ece B ^2 
q 



bi(C,q)+b a (C,-q) 



cosq y [ b^D, -q) + b a (D,q) ) + cosq^J D <-> B 
I 



, a = y, z. (25d) 



Here, the operators b(m, q) are the Fourier transforms of 
the bosonic operators b(m, r) 

b{m,q) ^e'' r i(m,r), m = A,B,C,D (26) 



where the summation of r is over the lattice of unit cells. 

This effective Hamiltonian has four decoupled parts, 
, H^ff , Hyff and H%" ■ Among these four decou- 
pled parts, Hyff and are just the Hamiltonian for 
the gapped magnetic mode on sublattice A and gapless 
magnetic Goldstone modes on sublattice B, C and D due 
to the spontaneous symmetry breaking of SO (3) spin 
symmetry. These two gapless Goldstone modes corre- 
spond to the gapless uniform rotation of the AF ordering 
from ^-direction to y or z direction. H^' is the Hamil- 
tonian of gapped x-magnon mode on sublattice B, C and 
D due to the condensation of x-magnons on these sub- 
lattices. The remaining Hf" is of great interest for it 
is the Hamiltonian of the charge modes and x-magnon 
mode on sublattice A. 

Under the insulating lobe of checkerboard state, the 
charge modes are gapped. The charge gap luq is given by 



the solution to the quartic equation of the form: 

w$ + C 2 ujI + c 3 = 0, 



(27) 



where C2 and C4 are quadratic functions of J c . This 
equation is quartic because fields with momentum q and 
— q are coupled and there are tH^ terms in the effective 
Hamiltonian. 

The condition C3 = gives the lobe-shaped second 
order phase boundary where one charge mode becomes 
soft and a phase transition from the PDW state to the 
SC state occurs. However, one must be careful since it 
is possible that the system takes a first-order transition 
before the charge mode softens. This indeed happens for 
lower left part of the lobe, where the PDW state takes 
a first order transition to a AF+SC mixed state. In the 
following discussion, we assume the tip and some of the 
left part of the lobe survive, as in the case plotted in Fig. 

m 

As one approaches the left part of the second order 
lobe following a trajectory with constant J c , the energy 
cost of removing hole-pair decreases and becomes zero 
at the phase boundary. This observation leads to the 
conclusion that a particle-like charge mode becomes soft 
on the left part of the lobe. Similarly, one can argue 
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FIG. 10: Dispersion of four gapless modes at a quantum 
multi-critical point, located tip of the 5 = 1/8 checkerboard 
lobe illustrated in Fig. |S| The two spin wave modes have ex- 
actly the same dispersion due to the SO (3) symmetry of AF 
moment. For the parameters shown in the figure, these four 
modes have the same speed at the multi-critical point. These 
four modes correspond the four Goldstone modes required by 
the quantum dynamical 50(5) symmetry. 



that a hole-like charge mode becomes soft at the right 
part. At the tip, where the left part and right part of the 
lobe meet, both particle-like and hole-like modes become 
soft. Consequently, an effective particle- hole symmetry 
is dynamically restored. One can check that O2 indeed 
vanishes at the tip of the lobe. 

At the tip of the lobe, the speed of gapless charge 
modes is determined by the interaction V c and V"/. For 
appropriate values of V c and V' c , the speed of the gap- 
less charge mode can be the same as the speed of the 
two gapless AF spin wave modes. When quantum 5*0(5) 
symmetry is spontaneously broken to a 50(4) symmetry, 
there should be exactly four degenerate gapless Gold- 
stone bosons. This model shows that such a dynamical 
symmetry is possible at the multi-critical point around 
the tip of the 5=1/8 lobe. Fig. ^1] plots the disper- 
sions of four gapless modes at the tip of AF insulating 
checkerboard lobe in the phase diagram of Fig. [5] 



V. PHASE DIAGRAM OBTAINED FROM 
QUANTUM MONTE-CARLO SIMULATION 

A. Numerical simulation 

Because of the bosonic nature, the minus-sign problem 
is absent in the quantum 50(5) model. Therefore, sim- 
ulations can be carried out for systems with sizes much 
larger than the ones available with fermionic QMC simu- 
lation in the physically interesting region. The pioneering 



numerical wor ks22i2ii22i2&Z4iS show that the projected 
50(5) model can give a realistic description of the phase 
diagram of the HTSC cuprates and account for many of 
their physical properties. In this section, we shall present 
the simulation of the 50(5) model with extended in- 
teractions using the Stochastic Series Expansion (SSE) 
method£Si2L2& with operator-loop update^. This Quan- 
tum Monte Carlo method was shown to be very efficient 
for simulations of harcore bosonic systema22iSiS. The 
overall topology of the phase diagram agrees well with 
the mean field calculation presented in the previous sec- 
tion, although the parameters are strongly renormalized. 

From the values given in section [n] we see that we can 
safely neglect which is rather small so that only re- 
mains J c and J s which are both positive. In this section, 

In order to avoid the notorious sign problem in the 
Quantum Monte Carlo simulations of the 5*0(5) model 
with extended interactions, all off-diagonal terms should 
be positive. On a square lattice with only nearest- 
neighbour non-zero off-diagonal terms, the sign of these 
matrix elements can be safely changed by a harmless uni- 
tary transformation acting on hopping terms in only one 
of the sublattices. 

For each simulation, the number of loops (or "worms" ) 
made during the loop operator updated is calculated self 
consistently during the thermalization part, such that on 
average the number of vertices visited by worms during 
each loop operator update is equal to C(n). Here (n) is 
the average number of non-Identity vertices in the oper- 
ator string (see Ref. [j^) and C a proportionnality con- 
stant, usually taken between 1 and 5 (the larger O, the 
more autocorrelations between successive Monte Carlo 
configurations are reduced). 

In order to plot the phase diagram, we should com- 
pute the order parameters corresponding to AF, SC and 
PDW phase. We use the superfluid density p c to locate SC 
phases. Indeed, p c can be related to the winding numbers 
of the world-lines which can be directly computed dur- 
ing the QMC simulations&SP.. Here the winding number 
only involves the charged particles, i.e. the hole pair hop- 
ping. We can take a similar definition with the magnetic 
bosons to define the spin stiffness. 

It is straightforward to measure the density-density 
correlations and their Fourier transform, the structure 
factors N(q x ,q y ) which indicate PDW phases : 



N(q) 



^y^^exp{iq.r){n h (i)n h (i 



r)) 



These quantities characterize the diagonal long-range or- 
der. On finite clusters, the structure factor at the appro- 
priate momentum diverges as the volume of the system 
in the ordered phase, so that, by plotting N(q x ,q y )/L 2 vs 
1/L {L is the linear size of the system), a scaling analysis 
can demonstrate long range order. 

Due to the intrinsic complexity of the projected SO(5) 
model and to the large number of interaction types, we 
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restricted the simulations to small lattices (4x4 up to 
10 x 10 lattices) at low temperatures (typically (3 = 10), 
to mainly be in the ground state. Even using the pow- 
erful SSE method, we found that for specific points in 
the phase diagram (near phase boundaries for example), 
autocorrelation times of different obscrvables or tunnel- 
ing times between two neighbouring phases can be long, 
decreasing the statistical precision. This sometimes pre- 
vents us from providing definitive statements about the 
nature of phase transitions. However, outside of these 
regions, we can clearly distinguish the nature of the dif- 
ferent phases. 



B. Limiting cases 

The pure projected SO (5) model corresponds to V c = 
V' c = V„ = and J c — 2J S . This model has been studied 
with the same technique^. A first-order transition be- 
tween AF and SC phases was observed. It was already 
recognized that a small V c and V c were enough to turn 
this transition into a second-order phase transition. 

Another interesting case occurs when the triplet den- 
sity becomes small. In that limit, the model reduces to 
one type of hard-core boson with hopping J c and nearest 
and next-nearest neighbor repulsion V c and V' c respec- 
tively. 



H = - J c^[4«^0') + #-' 



(ij) 



n h (i)n h {j) (28) 



so that we can fix J c — 1 as the energy unit. This model 
has been extensively studied in Ref. |63| and some results 
have been obtained with the same SSE method. Let us 
review a few useful results. 



1. Half-filling 

The phase diagram at half-filling is well known and 
exhibits solid phases with either (tt, 0) (stripe) or (tt, n) 
(checkerboard) structures. In between exists a superfluid 
phase with a non-zero superfluid density p c . We recover 
the same results as Hebert et alM : with our choice of 
interactions, by varying J c , we can drive the system from 
a superfluid toward a (tt, 0) (stripe) phase. 



2. Away from half-filling 

Away from half-filling, our grand-canonical algorithm 
is able to check whether we are in a phase separated state 
or not by looking at the histogram of the density during 
the simulation. On general grounds, the presence of a 
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FIG. 11: Fhe global phase diagram obtained by QMC. The 
parameters used in simulation are A s = 4.8, V c ~ 4.1010, 
V' c = 3.6329 and J n = Vx = 0. The lines are guides to the 
eye only. The overall topology of the phase diagram agrees 
well with the mean-field phase diagram of Fig. |SJ although 
the parameters are renormalized. 



double peak structure shows the presence of true phase 
separation in the system. 

With our choice of parameters, we find that away from 
the striped state, there is a phase, close to half-filling, 
where both N(0,ir)/ L 2 (PDW order parameter) and p c 
are finite, that is a supersolid phase. Moreover, there is 
no sign of phase separation so that it is a true homoge- 
neous phase, as claimed in Ref. |63[ . 



C. Global phase diagram 

Now that the parameters are fixed, we can compute 
the phase diagram for various J c and chemical potential 
A c (see Fig.ITT)). In order to discuss these results, we plot 
in the following most of the data as a function of doping 
5 = rifc/2 which depends on the chemical potential A c 
as shown on Fig. ^1 but has no strong finite-size effects. 
Let us comment on some results. 



1. Large J c 

For all values of J c , the superfluid density increases 
linearly with doping, for small doping, thus capturing a 
key piece of the Mott physics. For large J c (J c = 1.5 
for example), reducing hole density starting from S — 
1/4, we have a smooth decrease of superfluid density (sec 
Fig. ll3|) . At the same time, magnon density increases and 
gives rise to AF. Fig.lHlshows typical data of N(tt, 0)/L 2 
which is the PDW order parameter for various sizes. In 
order to get information on the thermodynamic limit, we 
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FIG. 12: S = rih/2 vs chemical potential for J c — 1.5 and 
J c = 0.6. The presence of plateaux indicates incompressible 
PDW phases. The finite size effects are rather small. 



FIG. 14: Scaling of (tt, 0) PDW order parameter (N(n,0)/L 2 ) 
for various sizes (L = 6 to 12 from top to bottom) with J c = 
1.5. We checked that it extrapolates to for all 6. 



2. Small J c 




0.00 0.05 0.10 0.15 0.20 0.25 



FIG. 13: Superfluid density as a function of doping for various 
J c . The "class Bl" trace with J c = 1.0 bears great similarity 
with the famous dome-like phase diagram of many cuprate 
families. Also, the "class B3"-like trace with J c < 0.8 has a 
pronounced dip at 8 = 1/8 as the phase diagram of LSCO 
family does. 



have performed a finite-size scaling of our data, and it 
vanishes for all fillings. We have a superconducting state 
for all dopings, with a possible coexistence with AF at 
low-doping. When the doping vanishes, we recover a pure 
AF phase. 



For smaller J c , the long-range interactions between 
hole pairs start to play a role. 

Large doping — The first example is the appearance 
of a (it, 0) PDW insulating phase at n h = 1/2 (8 = 1/4). 
Indeed, close to this doping, the magnon density is very 
small so that the model is similar to the single hardcore 
boson case. We therefore recover a transition at fixed 
ri/j = 1/2 (5 = 1/4) as a function of J c between a super- 
fluid phase (SC) and an insulating stripe phase (with a 
finite PDW order and a vanishing superfluid density). A 
finite-size scaling of the PDW order parameter showing 
a finite value is shown on Fig.^|for rih — 0.5 (6 = 0.25) 
and J c — 0.6. On Fig. El we see that below J c ~ 1.3, 
the superfluid density vanishes at nh — 1/2 (S = 1/4) 
in agreement with what has been found for the single- 
boson modelS. Mean-field results find a higher value 
(F c '/2~1.8). 

Intermediate doping — A second example of the inter- 
action effect is given by the appearance of an insulating 
PDW phase at n h = 0.25 (5 = 0.125) for low enough J c 
(Jc < 0.7). A finite PDW order parameter and p c = 
are shown again on Fig. 1151 and 1131 for J c = 0.6. This 
state corresponds to a localization of one hole-pair every 
4 sites, so that this = 1/4 (6 = 1/8) checkerboard 
also possesses a finite N(n,ir)/L 2 as shown on Fig. IT^l 
We find that this nt = 1/4 (S = 1/8) checkerboard is 
insulating, with a vanishing superfluid density. However, 
it could be possible to have a supersolid phase^i. 

For intermediate dopings 0.25 < < 0.5(0.125 < S < 
0.25) and small J c , we find a finite superfluid density 
and possibly a finite PDW order, that is a supersolid 
region. We checked that this phase is stable against phase 
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FIG. 15: (tt, 0) PDW order as a function of 8 for various 
sizes (L = 6, 8 and 10) for Jc = 0.6. We checked that, within 
the system sizes simulated, finite PDW order survives in the 
thermodynamic limit for n% > 0.25 (5 > 0.125). 



FIG. 16: (tt, 7t) PDW order as a function of 8 for various 
sizes (L = 4, 6, 8 and 10) and J c = 0.7. We checked that it is 
finite only for nn = 1/4 (8 = 1/8) checkerboard. 



separation. We cannot be conclusive about its extension 
in the phase diagram, due to the long autocorrelation 
time for some observables. This explains the large error 
bars in some parts of the phase diagram. However, since 
we know that such a phase exists and is stable close to 
half-filling, we can assume that it has a finite extension, 
ft would be interesting to check whether this supersolid 
also exists for a single hard-core boson model close to 
n,h = 1/4 (6 = 1/8), which is much easier to simulate. It 
is remarkable that for J c values close to the transition, 
we see a dip in the superfluid density fFig. [T^|> . which is 
due to the proximity to the insulating state. 
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3. Magnetic properties 

As shown on Fig. I17I the spin stiffness shows a mono- 
tonic decrease with chemical potential (or doping) so 
that it is difficult to locate precisely where it vanishes. 
For small doping 5 or chemical potential — A c , we ob- 
serve a rapid linear decrease so that we can estimate 
roughly where AF could vanish in the thermodynamic 
limit. These phase boundaries are shown on Fig. Illl and 
are in agreement with what has been found at the mean- 
field level. 

For small J c , AF seems to vanish at rih = 1/4 (S = 
1/8). However, our data of p s also show a shallow peak 
above this filling, which could indicate a reentrance of 
AF as found in mean-field. Unfortunately, with our cur- 
rent limitation on available sizes and statistics, we cannot 
conclude for sure about this possibility. 



FIG. 17: The spin stiffness versus chemical potential for dif- 
ferent sizes and J c . 

D. Summary 

The qualitative features found at the mean-field level 
are still present in a full numerical calculation and we 
have a very nice overall agreement (see Fig. Iffl and lllfl . 
Of course, exact critical values of J c for superconducting- 
insulator transition are different from mean-field values 
but this does not change the physical conclusions. 

It is very difficult to point precisely where AF vanishes 
since spin stiffness does not show any sharp drop. How- 
ever, we clearly see a strong reduction of p s , which seems 
to decrease linearly with doping. With the chosen param- 
eters, it seems plausible that both this line and the PDW 
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transition line merge close to the tip of the n = 1/4 lobe 
as was found at the mean-field level. This result was as- 
sociated to the dynamical restoration of SO (5) symmetry 
at this point. Our data might be an indication that this 
is indeed the case (we have taken the same V c and V' c so 
it seems pretty robust). However, a complete answer can 
only be provided by computing dynamical correlations, 
which is more involved and require good statistics. 



VI. EXPERIMENTAL CONSEQUENCE AND 
PREDICTIONS 

In this work we have constructed a single quantum 
Hamiltonian based on the projected SO(b) model with 
extended interactions, and presented detailed analytical 
and numerical calculations which give a consistent global 
phase diagram as depicted in Fig. 0] This schematic 
phase diagram is obtained from the quantitative model 
calculations summarized in Fig. El and 1111 This model 
captures the overall topology of the cuprates phase di- 
agram, including the dome-like feature of T c , which is 
determined within our model by the superfluid density, 
the 1/8 anomaly due to charge ordering, the coexistence 
of SC with AF and possibly with charge order. As men- 
tioned in the introduction, some of these features have 
been discussed in other theoretical contexts before, how- 
ever, it is rather remarkable that they are all reproduced 
by a single quantum model accessible by reliable QMC. 
Below we shall discuss some of these features and present 
more detailed theoretical predictions. 



A. Dependence of superfluid density on doping 

A remarkable feature of the HTSC cuprates is the 
dome-like dependence of T c on doping. Experiments 
have also shown a remarkable dependence of superfluid 
density on doping. On the underdoped side, both T c 
and the superfluid density scales linearly with doping, a 
fact commonly referred to as the "Uemura plot" . Fur- 
ther examination also shows that the superfluid density 
in the overdoped regime decreases with increasing dop- 
ing, which is commonly referred to as the "boomerang 
effect" . In cuprate system, the muon spin relaxation 
rate a(T — * 0) is proportional to n s /m* where n s is the 
superfluid density and m* is the effective mass of hole- 
pairs 82 . To explain the deviation from the linear relation 
between a and doping in the overdoped regime, Uemura 
proposed that some of the doped holes do not form pairs 
and are phase separated from the SC hole-pairs, even at 
zero temperature^. A similar phase separation picture 
was proposed by UchidaS. However, it has always been 
rather puzzling that if the phase separated normal elec- 
trons existed at zero temperature, they should provide 
a channel of relaxations, which has not been observed 
experimentally. 



Within our effective bosonic model, T c is directly de- 
termined by the superfluid density p c . As we can see 
from Fig. ED Pc indeed scales linearly with the doping 
density for small doping, and has a dome-like dependence 
for higher doping for intermediate values of J c . It peaks 
around 8 ~ 18% for J c = 1.2. The physical reason for the 
behavior on the overdoped side arises from the tendency 
of hole pairs to form a competing charge ordering state at 
either 3/16 or 1/4. Within this picture, when more holes 
are added into the system on the overdoped side, the 
holes are still paired, but some of them form a charge or- 
dered PDW state, with preferred doping of 3/16 or 1/4. 
This charge ordered PDW state either phase separate 
from the SC state, or coexists with it, but either way, 
the superfluid density is reduced because of the repul- 
sive coupling between the two forms of order. Therefore, 
this picture predicts a new charge ordered state on the 
overdoped side, which should be tested experimentally. 
Since the magnon density decreases monotonically with 
doping, the charge ordered states on the overdoped side 
may not be AF ordered, which makes it hard to observe 
by neutron scattering. Furthermore, on the overdoped 
side, our purely bosonic model also becomes less accu- 
rate, and a fully quantitative theory has to include the 
low energy fermionic excitations. 



B. The 1/8 anomaly and pressure experiments 

In the LSCO systems, T c has a pronounced dip near 
doping of S = 1/8 84 . More recently, it has been demon- 
strated that the competition between the nearly insulat- 
ing 1/8 phase and the SC phase in LSCO family can 
be controlled by pressure. There are two different ap- 
proaches. One is to apply hydrostatic or uniaxial pres- 
sure on single crystal^. In this way, the pressure in ab 
plane increases T c while the pressure in c direction de- 
creases T c . The other is to introduce compressive or ex- 
pansive strain into thin films with the help of the lattice 
mismatch between the film and substrate 7,1 . Enhanced 
T c and disappearance of the 1/8 anomaly for compressed 
film and strong reduction of T c for expanded film were 
observed^. 

Our model reproduces the 5=1/8 effect for J c < 1, 
as one can see from Fig. 1131 The pronounced dip is 
caused from the competition between the SC state and 
the insulating PDW states. Varying doping correspond 
to the "B2" or "B3" traces as depicted in Fig. 01 On the 
other hand, pressure in the ab plane reduces the lattice 
constant, thus increases the hopping term J c . Therefore, 
applying pressure in the ab plane is equivalent to follow- 
ing a "class A2" trace starting from a small J c in the 
global phase diagram. The doping dependence of T c for 
different films given in Fig. 3 of Sato et afii shares many 
common features with our doping dependence of SC order 
parameter, both MF results (Fig. [7J and QMC results 
(Fig. I13fl . The destructive effect on T c of the pressure 
in the c-direction can also be understood in terms of the 
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Poisson effect. The lattice constant in ab plane will in- 
crease when the sample is compressed in the c direction. 
This will lead to the decrease of T c as argued previously. 

Similar to the "class A2" trace, "class A3" trace can 
also be realized by applying pressure. Therefore, we 
predict a similar pressure induced superconductor-to- 
insulator transition at 5 = 1/16. 

In this way, the pressure effect on the 1/8 anomaly 
is understood in terms of the bosonic superfluid-to- 
insulator transition at the fixed density of 6 = 1/8. Stan- 
dard predictions on the superfluid-to-insulator transition 
applies to the 5 = 1/8 transition. In addition, we ar- 
gued that the tip of the 6=1/8 lobe can possibly have 
the full quantum 50(5) symmetry. This prediction can 
be experimentally tested by comparing both the static 
and dynamic charge and magnetic responses, as we have 
discussed in section II. B. 



compatible with the observed size of the charge unit cell. 
In Ref. 43$, the hole pair checkerboard state was estab- 
lished by a mean field calculation. A main result of the 
present work is to establish the existence of this state by 
QMC calculations. 

Our calculations as summarized in Fig. IB1 and HP show 
that the charge ordered insulating states are also accom- 
panied by AF magnetic order. Enhanced AF fluctuation 
in the vortex state was originally predicted within the 
SO (5) theory in Ref.[lllll7|. and has been experimentally 
observed in a number of HTSC cuprate families by a va- 
riety of experimental techniquesi£ii2i22i2ii22i2&2i. More 
recently, AF order has been detected by neutron scat- 
tering above H c2 , in electron doped cuprate a 25 ' 26 . This 
magnetic field induced quantum phase transition from 
the SC state to the AF state correspond to the "A2" 
trace as depicted in Fig. 0] 



C. The vortex phase and the ground state above 

The "class A2" and "class A3" traces can also be ap- 
proached in the cuprates by applying a magnetic field 
along the c-axis, which effectively reduces the hopping 
term J c . In underdoped LSCO samples, we predict that 
the magnetic field destroys SC order by localizing hole- 
pairs into a PDW state. This naturally leads to the field- 
induced insulating behavior in underdoped LSCO^^. 
For the YBCO and BSCO systems, the magnetic field 
could drive the hole-pairs into a disordered state before 
the lobe of 6 — 1/8 or 5 — 1/16 is reached. 

Recently, a striking feature is revealed in the STM ex- 
periments by Hoffman et alml, where the local density 
of states (DOS) near the vortex core show a two dimen- 
sional checkerboard- like modulation with a 4a x 4a charge 
unit cell. Here a is the lattice spacing of the C11O2 plane. 
This modulation decays exponentially away from the cen- 
ter of the vortex core, with a decay length of about 10a. 
A similar pattern has also been seen in the absence of the 
applied magnetic field^, possibly induced by the impuri- 
ties at the surface. In the case of vortex core experiment, 
SC is destroyed by the magnetic field, and the nature of 
the competing state is revealed. The case of impurity 
scattering is more complex, and the experimental obser- 
vation can also be interpreted as due to quasi-particle 
interference22i22i£i. 

The 6=1/8 insulating PDW state was proposed as an 
explanation for the STM measurements^. This state has 
the 4a x 4a checkerboard symmetry as observed in the ex- 
periment, and the doping level for the insulating 6=1/8 
state is reasonably close to the optimal doping level of 
the cuprates. On the other hand, if the holes themselves, 
rather than the hole pairs, form a Wigner crystal state, 
the periodicity of the charge ordering would be larger by 
a factor of y/2, inconsistent with the experiment. There- 
fore, by forming the Wigner crystal state of the hole pairs 
rather than the holes themselves, the doping level can be 



D. Charge localization and suppression of thermal 
conductivity 

In previous subsections, we argue that the competition 
between PDW and SC can be tuned by applying pres- 
sure or magnetic field. In particular, we show that the 
PDW ordering of hole-pairs can develop in the vortex 
core. The localized hole-pairs in a PDW state are ex- 
pected to have little contribution to the transport prop- 
erties such as thermal conductivity. This leads to the 
argument of the suppression of thermal conductivity by 
applying a magnetic field in c-axis below some tempera- 
ture Tk- Since the onset temperature Tk is expected to 
be proportional to the superfluid density, Tk has weak 
dependence on the magnetic field and follows T c in the 
underdoped cuprates. Moreover, the closer the system 
is to the charge ordered insulating state such as the 1/8 
PDW state, the smaller the suppression effect would be. 
Finally, the in-plane magnetic field has little effect on 
the thermal conductivity due to the fact that it does not 
create vortex in the ab plane. 

Therefore, trace "A2" in our global phase diagram and 
the charge localization into a hole pair crystal can pos- 
sibly explain the recent experiment on the suppression 
of thermal conductivity by applying a magnetic field in 
c-directionS. We also predict that applying pressure will 
also induce the suppression or enhancement of the ther- 
mal conductivity around the 1/8 doping, assuming the 
pressure will not induce strong effect of the lattice struc- 
ture which can also change the thermal conductivity. 

E. 1/16 doping 

As discussed in section ITTT1 additional insulating lobes 
at 5 = 1/16 and 6 = 3/16 doping levels are predicted 
if interactions are more extended than the nnn interac- 
tions included in this work. The charge and spin or- 
dering patterns for this state are depicted in Fig. ED 
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FIG. 18: The hole-pair checkerboard state with anti-phase 
magnetic domain at doping level S = 1/16. Circles denote 
holes and arrows denote spins while black crosses denote spins 
in a singlet bond. 
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FIG. 19: The hole-pair checkerboard state with in-phase mag- 
netic domain at doping level 5 = 1/8. Circles denote holes 
and arrows denote spins while black crosses denote spins in a 
singlet bond. 



FIG. 20: The hole-pair checkerboard state with anti-phase 
magnetic domain at doping level 8 = 1/8. Circles denote holes 
and arrows denote spins while black crosses denote spins in a 
singlet bond. 

netic order in the NCCO materialise, but not consis- 
tent with the magnetic structure observed in the LSCO 
material, which has anti-phase domain walls. We would 
like to stress that this is not a limitation on the fun- 
damental approach taken here. As shown in Fig. [201 
checkerboard charge order of the hole pairs can be fully 
consistent with the AF magnetic structure with anti- 
phase domain walls. However, stability of this type of 
magnetic structure requires more extended magnetic in- 
teractions. Since the complexity of both the CORE algo- 
rithm and the QMC increases substantially, we have not 
yet been able to derive such extended interactions from 
the microscopic models, and simulate them with QMC. 



G. Coexistence phases 



Even though the current paper does not investigate this 
type of more extended interactions explicitly due to nu- 
merical complexities, it is clear that the physics of these 
PDW states are similar to the 1/8 doping. Preliminary 
evidence for the S = 1/16 insulating state exists for 
the LSCO materialman. As we see from Fig. [H and 
1191 the charge ordering pattern rotates from diagonal at 
S = 1/16 to square at 6 = 1/8. This transition between 
states with different rotational symmetries could possibly 
be related to the diagonal to square transition observed 
in the neutron scattering experiments at S ~ 5% in the 
LSCO material^. 



F. Magnetic order 

In this paper, we focused our discussion on the checker- 
board charge order at 8 = 1 /8. Strictly within our model, 
the accompanied AF magnetic order is commensurate, as 
sketched in Fig. ^5] This type of magnetic structure is 
consistent with the recently observed field induced mag- 



While the simple SO(5) modeU> predicts the coexis- 
tence phases of AF and SC, more regions with coexisting 
charge, AF and SC orders are predicted in the global 
phase diagram of Fig. and Fig. ^2 I n particular, the 
"type 1" first order transition from undoped AF state to 
SC state predicted by simple SO(5) model is turned into 
two "type 2" second order transitions by the interactions 
V c and V£, which is consistent with previous study^i. 
While the width of these coexistence regions is model 
dependent, coexistence phases are important qualitative 
predictions of our theory. Experiments on cuprates have 
indeed suggested such coexistence phase a 89 i 90 ' 91 i 92 ' 93 i 94 . 

As we see from Fig. |H] the AF order disappears around 
S = 0.10. This is consistent with the value obtained from 
the t-J modei 9 ^. The width of the AF/SC coexistence 
phase largely depends on the values of J w and V„. As 
determined by CORE method, J w is negative and V v 
is positive. Thus, at the mean-field level, both and 
Vn terms induce a repulsion between AF and SC in the 
mixed state of coexisting AF and SC. When these terms 
are included, one would expect a smaller region of the 
AF/SC coexistence state of in the global phase diagram. 
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On the other hand, these two terms have different effects 
on the checkerboard state. Since the mean-field value of 
7r-operator in the PDW state is zero, the interaction 
will not induce any interaction terms at the mean-field 
level. In contrast, the term effectively changes the 
local chemical potential of hole-pairs (magnons) due to 
the nonzero mean-field value of magnons (hole-pairs) on 
nn plaquette. Therefore, also reduces the height of 
the insulating PDW lobe. 

VII. CONCLUSIONS 

Starting from commonly used microscopic models of 
high T c cuprates, an effective bosonic model can be de- 
rived by the CORE algorithmS&SL. In addition to the 
simple interactions included in the original projected 
50(5) model2L2£Liiii, extended interactions play an im- 
portant role in determining the global phase diagram of 
the model. This model can be studied systematically 
by the analytical mean field theory and by the QMC 
method, thanks to the absence of the minus sign problem. 
The global phase diagram consists of self-similar insulat- 
ing PDW phases at rational filling fractions, immersed 



1 S. C. Zhang, Science 275, 1089 (1997). 

2 E. Dernier and S. C. Zhang, Phys. Rev. Lett. 75, 4126 
(1995). 

3 E. Dernier and S. C. Zhang, Nature 396, 733 (1998). 

4 J. Rossatmignod, L. P. Regnault, C. Vettier, P. Bourges, 
P. Burlet, J. Bossy, J. Y. Henry, and G. Lapertot, Physica 
C 185, 86 (1991). 

5 H. A. Mook, M. Yethiraj, G. Aeppli, T. E. Mason, and 
T. Armstrong, Phys. Rev. Lett. 70, 3490 (1993). 

6 H. F. Fong, B. Keimer, P. W. Anderson, D. Reznik, F. Do- 
gan, and I. A. Aksay, Phys. Rev. Lett. 75, 316 (1995). 

7 P. Dai, M. Yethiraj, H. A. Mook, T. B. Lindemer, and 

F. Dogan, Phys. Rev. Lett. 77, 5425 (1996). 

8 H. F. Fong, B. Keimer, D. Reznik, D. L. Milius, and I. A. 
Aksay, Phys. Rev. B 54, 6708 (1996). 

9 H. A. Mook, P. C. Dai, S. M. Hayden, G. Aeppli, T. G. 
Perring, and F. Dogan, Nature 395, 580 (1998). 

10 P. C. Dai, H. A. Mook, and F. Dogan, Phys. Rev. Lett. 
80, 1738 (1998). 

11 H. F. Fong, P. Bourges, Y. Sidis, L. P. Regnault, J. Bossy, 

A. Ivanov, D. L. Milius, I. A. Aksay, and B. Keimer, Phys. 
Rev. B 61, 14773 (2000). 

12 H. F. Fong, P. Bourges, Y. Sidis, L. P. Regnault, A. Ivanov, 

G. D. Gu, N. Koshizuka, and B. Keimer, Nature 398, 588 
(1999). 

13 H. He, Y. Sidis, P. Bourges, G. D. Gu, A. Ivanov, 
N. Koshizuka, B. Liang, C. T. Lin, L. P. Regnault, 
E. Schoenherr, et al., Phys. Rev. Lett. 86, 1610 (2001). 

14 H. He, P. Bourges, Y. Sidis, C. Ulrich, L. P. Regnault, 
S. Pailhes, N. S. Berzigiarova, N. N. Kolesnikov, and 

B. Keimer, Science 295, 1045 (2002). 

15 C. Stock, W. J. L. Buyers, R. Liang, D. Peets, Z. Tun, 
D. Bonn, W. H. Hardy, and R. J. Birgeneau, Phys. Rev. 



in the background of the uniform SC phase, as depicted 
schematically in Fig. Different families of cuprates are 
attributed to different traces in the global phase diagram. 
The overall topology of the phase diagram obtained from 
our model agrees well with the experiments, and so are 
the behaviors of various physical quantities. Inclusion of 
longer ranged interactions could bring detailed and quan- 
titative agreements with the cuprate phase diagram. 



VIII. ACKNOWLEDGMENTS 

We would like to acknowledge useful discussions with 
Drs. A. Auerbach, E. Demler, W. Hanke, J. P. Hu, 
S. Kivelson, E. Mukamel, G. Schmid, M. Troyer and 
C.J. Wu. This work is supported by the NSF under 
grant numbers DMR-9814289, and the US Department 
of Energy, Office of Basic Energy Sciences under contract 
DE-AC03-76SF00515 and by the Swiss National Science 
Foundation. HDC is also supported by a Stanford Grad- 
uate Fellowship. Part of the simulations were performed 
on the Asgard beowulf cluster at ETH Zurich and at 
SLAC. 



B 69, 014502 (2004). 

16 J. M. Tranquada, C. H. Lee, K. Yamada, Y. S. Lee, L. P. 
Regnault, and H. W. Ronnow, Phys. Rev. B 69, 174507 
(2004). 

17 D. P. Arovas, A. J. Berlinsky, C. Kallin, and S. C. Zhang, 
Phys. Rev. Lett. 79, 2871 (1997). 

18 S. Katano, M. Sato, K. Yamada, T. Suzuki, and T. Fukase, 
Phys. Rev. B 62, 14677 (2000). 

19 V. F. Mitrovic, E. E. Sigmund, M. Eschrig, H. N. Bach- 
man, W. P. Halperin, A. P. Reyes, P. Kuhns, and W. G. 
Moulton, Nature 413, 501 (2001). 

20 B. Lake, G. Aeppli, K. N. Clausen, D. F. McMorrow, 
K. Lefmann, N. E. Hussey, N. Mangkorntong, M. Nohara, 
H. Takagi, T. E. Mason, et al., Science 291, 1759 (2001). 

21 B. Lake, H. M. Ronnow, N. B. Christensen, G. Aep- 
pli, K. Lefmann, D. F. McMorrow, P. Vorderwisch, 
P. Smeibidl, N. Mangkorntong, T. Sasagawa, et al., Na- 
ture 415, 299 (2002). 

22 R. I. Miller, R. F. Kieft, J. H. Brewer, J. E. Sonier, 
J. Chakhalian, S. Dunsiger, G. D. Morris, A. N. Price, 
D. A. Bonn, W. H. Hardy, et al., Phys. Rev. Lett. 88, 
137002 (2002). 

23 B. Khaykovich, Y. S. Lee, R. W. Erwin, S. H. Lee, S. Waki- 
moto, K. J. Thomas, M. A. Kastner, and R. J. Birgeneau, 
Phys. Rev. B 66, 014528 (2002). 

24 V. F. Mitrovic, E. E. Sigmund, W. P. Halperin, A. P. 
Reyes, P. Kuhns, and W. G. Moulton, Phys. Rev. B 67, 
220503 (2003). 

25 H. J. Kang, P. C. Dai, J. W. Lynn, M. Matsuura, J. R. 
Thompson, S. C. Zhang, D. N. Argyriou, Y. Onose, and 
Y. Tokura, Nature 423, 522 (2003). 

26 M. Fujita, A. Matsuda, S. Katano, and K. Yamada, cond- 
mat/0311269 (unpublished). 



19 



27 S. C. Zhang, J. P. Hu, E. Arrigoni, W. Hanke, and A. Auer- 
bach, Phys. Rev. B 60, 13070 (1999). 

28 E. Altman and A. Auerbach, Phys. Rev. B 65, 104508 
(2002). 

29 S. Capponi and D. Poilblanc, Phys. Rev. B 66, 180503 
(2002). 

30 A. Dorneich, W. Hanke, E. Arrigoni, M. Troyer, and S. C. 
Zhang, Phys. Rev. Lett. 88, 057003 (2002). 

31 J. A. Riera, Phys. Rev. B 66, 134523 (2002). 

32 J. A. Riera, Phys. Rev. B 65, 174526 (2002). 

33 M. Jostingmeier, E. Arrigoni, W. Hanke, and S.-C. Zhang, 
Phys. Rev. B 68, 245111 (2003). 

34 J. M. Tranquada, B. J. Sternlieb, J. D. Axe, Y. Nakamura, 
and S. Uchida, Nature 375, 561 (1995). 

35 G. Aeppli, T. E. Mason, S. M. Hayden, H. A. Mook, and 
J. Kulda, Science 278, 1432 (1997). 

36 B. O. Wells, Y. S. Lee, M. A. Kastner, R. J. Christianson, 
R. J. Birgeneau, K. Yamada, Y. Endoh, and G. Shirane, 
Science 277, 1067 (1997). 

37 J. E. Hoffman, E. W. Hudson, K. M. Lang, V. Madhavan, 
H. Eisaki, S. Uchida, and J. C. Davis, Science 295, 466 

(2002) . 

38 C. Howald, H. Eisaki, N. Kaneko, M. Greven, and A. Ka- 
pitulnik, Phys. Rev. B 67, 014533 (2003). 

39 K. McElroy, R. W. Simmonds, J. E. Hoffman, D. H. Lee, 
J. Orenstein, H. Eisaki, S. Uchida, and J. C. Davis, Nature 
422, 592 (2003). 

40 Q. H. Wang and D. H. Lee, Phys. Rev. B 67, 020511 

(2003) . 

41 K. McElroy, J. E. Hoffman, D. H. Lee, K. M. Lang, 
H. Eisaki, S. Uchida, and J. C. Davis, Physica C 388, 
225 (2003). 

42 S. A. Kivelson, I. P. Bindloss, E. Fradkin, V. Oganesyan, 
J. M. Tranquada, A. Kapitulnik, and C. Howald, Rev. 
Mod. Phys. 75, 1201 (2003). 

43 H. D. Chen, J. P. Hu, S. Capponi, E. Arrigoni, and S. C. 
Zhang, Phys. Rev. Lett. 89, 137004 (2002). 

44 S. Kivelson, D. H. Lee, and S. C. Zhang, Phys. Rev. B 46, 
2223 (1992). 

45 L. P. Pryadko and S. C. Zhang, Phys. Rev. B 54, 4953 
(1996). 

46 C. P. Burgess and B. P. Dolan, Phys. Rev. B 63, 155309 
(2001). 

47 E. Witten, hep-th/0307041 (unpublished). 

48 H. Chen, O. Vafek, A. Yazdani, and S. Zhang, cond- 
mat/0402323 (unpublished). 

49 J. Zaanen and O. Gunnarsson, Phys. Rev. B 40, 7391 
(1989). 

50 H. Tsunetsugu, M. Troyer, and T. M. Rice, Phys. Rev. B 
51, 16456 (1995). 

51 S. A. Kivelson, E. Fradkin, and V. J. Emery, Nature 393, 
550 (1998). 

52 S. R. White and D. J. Scalapino, Phys. Rev. Lett. 80, 1272 
(1998). 

53 J. Zaanen, Science 286, 251 (1999). 

54 V. J. Emery and S. A. Kivelson, Nature 374, 434 (1995). 

55 M. Vojta and S. Sachdev, Phys. Rev. Lett. 83, 3916 (1999). 

56 E. Dernier, S. Sachdev, and Y. Zhang, Phys. Rev. Lett. 87, 
067202 (2001). 

57 Y. Zhang, E. Dernier, and S. Sachdev, Phys. Rev. B 66, 
094501 (2002). 

58 R. B. Laughlin, cond-mat/0209269 (unpublished). 

59 S. Haas, E. Dagotto, A. Nazarenko, and J. Riera, Phys. 
Rev. B 51, 5989 (1995). 



60 Y. H. Kim and P. H. Hor, Mod. Phys. Lett. B 15, 497 

(2001) . 

61 M. P. A. Fisher, P. B. Weichman, G. Grinstein, and D. S. 
Fisher, Phys. Rev. B 40, 546 (1989). 

62 C. Pich and E. Frey, Phys. Rev. B 57, 13712 (1998). 

63 F. Hebert, G. G. Batrouni, R. T. Scalettar, G. Schmid, 
M. Troyer, and A. Dorneich, Phys. Rev. B 65, 014513 

(2002) . 

64 C. J. Morningstar and M. Weinstein, Phys. Rev. D 54, 
^ 4131 (1996). 

65 S. Chakravarty, B. I. Halperin, and D. R. Nelson, Phys. 
Rev. Lett. 60, 1057 (1988). 

66 E. Manousakis, Rev. Mod. Phys. 63, 1 (1991). 

67 A. Vanotterlo and K. H. Wagenblast, Physical Review Let- 
ters 72, 3598 (1994). 

68 G. G. Batrouni, R. T. Scalettar, G. T. Zimanyi, and A. P. 
Kampf, Physical Review Letters 74, 2527 (1995). 

69 J. M. Kosterlitz, D. R. Nelson, and M. E. Fisher, Phys. 
Rev. B 13, 412 (1976). 

70 S. Arumugam, N. Mori, N. Takeshita, H. Takashima, 
T. Noda, H. Eisaki, and S. Uchida, Phys. Rev. Lett. 88, 
247001 (2002). 

71 H. Sato, A. Tsukada, M. Naito, and A. Matsuda, Phys. 
Rev. B 62, R799 (2000). 

72 D. G. Hawthorn, R. W. Hill, C. Proust, F. Ronning, 
M. Sutherland, E. Boaknin, C. Lupien, M. A. Tanatar, 
J. Paglione, S. Wakimoto, et al., Phys. Rev. Lett. 90, 
197004 (2003). 

73 X. F. Sun, S. Komiya, J. Takeya, and Y. Ando, Phys. Rev. 
Lett. 90, 117004 (2003). 

74 X. Hu, J. Low Temp. Phys. 117, 289 (1999). 

75 X. Hu, Phys. Rev. Lett. 87, 057004 (2001). 

76 A. W. Sandvik and J. Kurkijarvi, Phys. Rev. B 43, 5950 
(1991). 

77 A. W. Sandvik, Phys. Rev. B 56, 11678 (1997). 

78 A. W. Sandvik, Phys. Rev. B 59, 14157 (1999). 

79 G. Schmid, S. Todo, M. Troyer, and A. Dorneich, Phys. 
Rev. Lett. 88, 167208 (2002). 

80 E. L. Pollock and D. M. Ceperley, Phys. Rev. B 36, 8343 
(1987). 

M. Troyer, private communication. 

82 Y. J. Uemura, Solid State Commun. 120, 347 (2001). 

83 S. Uchida, Solid State Commun. 126, 57 (2003). 

84 A. R. Moodenbaugh, Y. W. Xu, M. Suenaga, T. J. Folkerts, 
and R. N. Shelton, Phys. Rev. B 38, 4596 (1988). 

85 K. Kudo, M. Yamazaki, T. Kawamata, T. Adachi, 
T. Noji, Y. Koike, T. Nishizaki, and N. Kobayashi, cond- 
mat/0311402 (unpublished). 

86 F. Zhou, P. H. Hor, X. L. Dong, W. X. Ti, J. W. Xiong, 
and Z. X. Zhao, cond-mat/0309034 (unpublished). 

87 M. Fujita, K. Yamada, H. Hiraka, P. M. Gehring, S. H. Lee, 
S. Wakimoto, and G. Shirane, Phys. Rev. B 65, 064505 
(2002). 

88 M. Matsuda, M. Fujita, K. Yamada, R. J. Birgeneau, M. A. 
Kastner, H. Hiraka, Y. Endoh, S. Wakimoto, and G. Shi- 
rane, Phys. Rev. B 62, 9148 (2000). 

89 J. H. Brewer, J. F. Carolan, W. N. Hardy, B. X. Yang, 
P. Schleger, R. Kadono, J. R. Kempton, R. F. Kiefl, S. R. 
Kreitzman, G. M. Luke, et al., Physica C 162, 33 (1989). 

90 C. Niedermayer, C. Bernhard, T. Blasius, A. Golnik, 
A. Moodenbaugh, and J. I. Budnick, Phys. Rev. Lett. 80, 
3843 (1998). 

91 J. E. Sonier, J. H. Brewer, R. F. Kiefl, R. I. Miller, G. D. 
Morris, C. E. Stronach, J. S. Gardner, S. R. Dunsiger, 



20 



D. A. Bonn, W. N. Hardy, et al., Science 292, 1692 (2001). 
Y. Sidis, C. Ulrich, P. Bourges, C. Bernhard, C. Nieder- 
mayer, L. P. Regnault, N. H. Andersen, and B. Keimer, 
Phys. Rev. Lett. 86, 4100 (2001). 

H. A. Mook, P. C. Dai, S. M. Hayden, A. Hiess, J. W. 
Lynn, S. H. Lee, and F. Dogan, Phys. Rev. B 66, 144513 



(2002). 

R. I. Miller, R. F. Kiefl, J. H. Brewer, J. E. Sonier, D. A. 
Bonn, W. H. Hardy, R. Liang, J. M. Mnard, and P. Poon, 
unpublished. 

A. Himeda and M. Ogata, Phys. Rev. B 60, R9935 (1999). 



